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Abstract 

In this paper, we provide the different types of bifurcation diagrams for a super- 
conducting cylinder placed in a magnetic field along the direction of the axis of the 
cylinder. The computation is based on the numerical solutions of the Ginzburg-Landau 
model by the finite element method. The response of the material depends on the val- 
ues of the exterior field, the Ginzburg-Landau parameter and the size of the domain. 
The solution branches in the different regions of the bifurcation diagrams are analyzed 
and open mathematical problems are mentioned. 

1 Introduction 

In this paper, we want to give a detailed description of the bifurcation diagrams of an 
infinite superconducting cylinder of cross section submitted to an exterior magnetic 
field ho. The response of the material varies greatly according to the value of ho, the 
size of the cross section and the Ginzburg-Landau parameter k that characterizes the 
material: superconductivity appears in the volume of the sample for low fields and small 
samples, under the form of vortices for higher fields, bigger samples and larger values of 
k, and is destroyed for high fields. The type of response of a superconducting material 
has been studied numerically and theoretically by various authors in various asymptotic 

regimes §, |o|, n, m I]], m, m n> m m • Here ' 

we want to give a complete numerical 
picture of the bifurcation diagrams for all values of the parameters. This type of 
computation has already been done in dimension 1 by one of the authors using auto 
P|. Here, the computation is made in a 2 dimensional domain using numerical solutions 



of the well-known Ginzburg-Landau model [31], based on a code first developed in Jl 



We examine the behavior of the energy, the magnitude of the order parameter and the 
magnetization versus the magnetic field for various solution branches. We also provide 
some analysis on the detailed findings such as the stability of solutions. 
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The paper is organized as follows: the Ginzburg-Landau model is briefly stated in 
Section 2 and the main features of the numerical codes used in the computation are 
described Section 3. The complete phase diagrams are given in Section 4, along with 
detailed analysis. A conclusion is given in Section 5. 



2 The Ginzburg-Landau model 

Let f2 denote a two dimensional bounded domain which represents the cross section of 
a three dimensional cylinder occupied by the superconducting sample. Assume that 
the cylinder is homogeneous along its axis and a constant applied field Hq is placed 
along the axis direction as well. Then, the Gibbs free energy Q may be written in the 
following form J31[|: 



dn. 



S(V^)=^(/ n + a|Vf + fM 4 ) dn 



2iu f, 



|curl A\ 2 _ curl A-H 



Here, ip is the (complex-valued) order parameter, A is the magnetic potential, curl A is 
the magnetic field, a and (3 are constants (with respect to the space variable x) whose 
values depend on the temperature, c is the speed of light, e s and m s are the charge 
and mass, respectively, of the superconducting charge-carriers, and 2^% is Planck's 
constant. 

After proper nondimensionalization, we can reformulate the free energy functional 

as: 

Q(ip,A) = J |(V - iA) 4>\ 2 + y (l - M 2 ) 2 + Icurl A - h \ 2 dn, 

where k is the Ginzburg-Landau parameter representing the ratio of the penetration 
depth and the coherence length, ho the applied magnetic field and d the characteristic 
size of the domain 0, that is = dD where D is a fixed domain. 

The system that we are going to study is the following Ginzburg-Landau equations 



derived as the Euler-Lagrange equations for the minimizers of the functional Q [21]: 

r _( V -a)> = K 2 v>(i-ivf) fcfi, m 

\ -curl curl A = \^\ 2 A + §(V>*VV> - i^V^*) in Q, { ' 

which are supplemented by the boundary conditions 

f (Vip - iAip) • n = on dtt, 
| curl A = /iq on dQ, 



and gauge constraints 



div A = in ft, 
A • n = on dQ. 



Here, 50 is the boundary of Q and n is its unit outer normal. With the above nondi- 
mensionalization, takes values between and 1: the normal state corresponds to 
|^| = while the Messiner state corresponds to = 1. 

This Ginzburg-Landau model has a special family of solutions called the normal 
solutions: ip = and curl A = /io, which correspond to the situation where supercon- 
ductivity is destroyed. According to the values of the different parameters k, d and 
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the system may have other solutions: superconducting solutions, for which if; is never 
and vortex solutions for which ip has isolated zeroes (see Figure |] and ||). For a 
complete introduction to the topic, one may refer to |31| |. 



3 Numerical codes 



The numerical code that we are using is based on the finite element approximations of 



the Ginzburg-Landau model, first proposed in [17] and later used in many settings, see 



1 18, 32, 23] for instance. The codes have many different variants that can be used to 
simulate the solutions of the Ginzburg-Landau models in the high-K setting flu], 16 1, 
thin films with variable thickness pi, samples with normal inclusions |18j as well as 



layered materials. Here, we choose the standard version that solves the Ginzburg- 
Landau equations on a rectangular domain f2 and d will denote the characteristic size 
of the rectangle. We note that other numerical methods, such as the gauge invari- 
ant difference methods fjl!|, i.e, the so-called bond-and-link variable methods or the 
method of eigen-functions, have also been used to compute the phase diagrams for the 
Ginzburg-Landau equations [||, 14 1. 

In our implementation of the finite element approximation, we use a uniform tri- 
angular grid with piecewise quadratic polynomials for both ip and A. It is shown that 



(see [17 1 for instance), if (ips,As) is the finite element solution on a given mesh with 
mesh size 5, the convergence of the approximation is assured, and asymptotically, we 
have 

\ 2 + \\A-A s \\ 2 = 0(5 3 ) 



where || • [[2 denote the standard mean square I? norm. For each set of calculation, we 
refine the mesh size until numerical convergence is evident. 

In Figures [l] and ^, we present a few typical plots for the numerical solutions of the 
equation (|l|). For each solution, the plots include a surface plot of the magnitude of the 
order parameter, a surface plot of the magnetic field given by curl A and a vector plot 
of the superconducting current. In Figure [l], we have a solution with a single vortex at 
the center of the domain which corresponds to the parameter values k = 0.23, d = 16.8 
and ho = 0.563. ... 





Figure 1: A single vortex solution of ([!]). 



In Figure we present the plots for a solution with two vortices corresponding to 
k = 0.8, d = 4 and h = 1.2. 
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Figure 2: A solution of (|l|) with two vortices. 



For fixed k, d and we are interested in finding the number of solutions of ((T|) and 
their stability. A continuation in the parameter spaces is used for getting solutions with 
different parameter values. With k, d given but ho allowed to vary, for a computed 
solution branch, we plot HV'Hoo, the maximum magnitude of the order parameter ip 
in the domain, the free energy Q and the magnetization versus the applied field Hq. 
These phase diagrams or bifurcation diagrams will give us information on the solutions 
(number and stability) for each h^. They were drawn by Ginzburg |2(J in some limiting 
cases of the parameters. Here we want to give a more complete description of these 
diagrams for all values of the parameters. 

The results of our numerical computations allow us to separate the n-d plane into 
different regions depending on the shape of the bifurcation diagram. A solution branch 
may exist in one region but may cease to exist in another one, it may also develop 
hysteresis in some regions. 



4 The bifurcation diagrams 

It is well known that for large fields {h > h*), the only solution is the normal solution 



1 22]. For smaller fields (h < h*), the normal solution always exists but there are other 
solutions which display four different types of behaviors. These behaviors depend on 
the values of k and d. In Figure ||, we have plotted four curve segments {Ki(d)}f 
separating the k — d plane into four regions {Ri}f- 

The four curve segments in Figure ^ are the results of the computation described 
in the earlier section. All four curves meet close to k = 1/V2, K2(d) is of the form 
2.112/d, Ks(d) is tending to 0.4 at infinity. The four regions correspond to the four 
types of behaviors for the bifurcation diagrams. For convenience, for each i = 1, 2, 3, 4, 
we use d = di{K) to denote the inverse function of the function k = Hi(d) wherever the 
inverse is well-defined. What distinguishes the different regions are features like the 
existence (or the lack of existence) of vortex solutions, the global and local stability of 
solutions, and the hysteresis phenomena. 

We now provide detailed descriptions of the solution behavior for each region in 
Figure [| 

Region 1: d < di{n) and d < ^(ft). 

This corresponds to the situation where the cross section of the superconducting sample 
is small enough. The bifurcation diagram is illustrated in Figure ||. The corresponding 
plot of the energy is given in Figure [5] and the magnetization curve in Figure ||. 
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Figure 3: The curves Ki(d), K2(d), n 3 (d) and K^{d) 
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Figure 4: The bifurcation curve for d = 2.0 k = 0.3 

Throughout this region, there is a unique non- normal solution for h < h*. This 
solution is a superconducting solution which is the global minimizer of the free energy 
Q. The curve l^loo against h is monotonically decreasing. When increasing the field, the 
magnitude of the superconducting solution decreases until it turns normal at ho = h*. 
Conversely, when decreasing the field, the normal solution turns superconducting also 
for ho = h*. The transition to the normal solution is of second order, that is the 
energy of the superconducting solution tends to the energy of the normal solution at 
the transition and there is no hysteresis phenomenon. This can also be seen on the 
magnetization curve of Figure |(| 

There is no vortex solution for the parameters (d, k) in this region. This reflects 
the fact that d is too small to allow enough room for a vortex to exist since a vortex 
core is of typical size C / n. 

In ||, it is proved that for d < minjdo, di/n}, a solution of the Ginzburg-Landau 
equations has no vortex, and for such a solution, ip is almost constant in the domain 
and is unique up to multiplication by a constant of modulus 1. Also it is proved that 
the bifurcation curve for |V>|oo is decreasing and that h* is asymptotically C/d when 
d tends to 0. More precisely, for a disc, when d tends to 0, the limiting curve is 
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Figure 5: The energy for d = 2.0 k = 0.3 
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Figure 6: The magnetization for d = 2.0 k = 0.3 

11^11^ = 1-^/2^. 

Region 2: k > K2{d) and k > ^(d) . 

In contrast with Region 1, this region corresponds to the situation where the typical 
size of vortices (C/k) is small enough compared to the size of the domain. This region 
displays the typical type II behavior of superconductors. The bifurcation diagram is 
illustrated in Figure 0, the energy in Figure |8| and the magnetization curve is given in 
Figure []. 

It has been well understood both physically and in more recent years mathemat- 
ically that, for sufficiently large k, there are vortex solutions which are the global 
minimizers of the free energy for a certain range of fields. The number of vortices de- 
pends on the strength of the applied field. The maximum number of vortices increases 
with d and type-II behavior. 

For very small fields, the global minimizer is the superconducting solution (solid 
line). As the field is increased, the superconducting solution loses its global stability 
(h = h Cl and for even larger fields loses its local stability, h = /iq). Then the global 
minimizer starts to nucleate vortices. In Figures 0, |, |, the solut ion branch corre- 



6 




sponding to a one vortex solution is illustrated by a dotted line and to a two vortex 
solution by a dashed line. 

Let n be a non- negative integer (the flux quanta). Hysteresis phenomena occur in 
our experiments when going from a solution with n vortices to a solution with n + 1 
vortices: when increasing the field, the solution with n vortices remains locally stable, 
though it is no longer a global minimizer until h* n where it jumps to a solution with n+1 
vortices. Similarly when decreasing the field, the solution with n + 1 vortices remains 
locally stable down to < h* n . Such a hysteresis pattern is typical for parameters 
in this region. The fields /i* and correspond to the critical field where the vortex 
solution loses its local stability and their values can be identified from Figures [?] and |8| 
(transition from superconducting solution to 1 vortex solution and from one vortex to 
two vortices). 

Nevertheless, the transition from normal to the last vortex solution is of second order 
since the vortex solution is stable, as can be seen from the energy and magnetization 
curves. This transition occurs at a field that is usually called h C3 in the literature. 

There have been many numerical works on these vortex solutions: see for instance 
1 13]. Recently, a rigorous analysis of vortices, their number and the critical thermody- 
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Figure 9: The magnetization for d = 3.2 k = 1 



namic fields in the high kappa limit has been done by Sandier and Serfaty p7| , 28, 29 1: 
in particular they obtain in the high kappa limit that h Cl (the field for which the energy 
of the superconducting solution is equal to the energy of the one vortex solution) is of 
order C log k and the field h* n is asymptotic to h Cl + n log log k. 

The onset of superconductivity in decreasing fields (instability of normal solutions 
and computation of the fields of nucleation) has been analyzed by Bernoff and Sternberg 



|8[ and del Pino, Felmer, Sternberg [12|. Their works provide, as d and n tends to oo, 
an asymptotic development of /i C 3, the field at which the normal solution bifurcates to 
a vortex solution. This is what is called surface superconductivity. A linearizing the 
Ginzburg-Landau equation has been done near the normal solution. In the high kappa 
limit, their computation yields 

K 2 

h C3 ~ — + CKK max + o(k) (2) 

Al 

where n max is the maximal curvature of the domain and Ai is the first eigenvalue of 
the linearized problem and is approximately equal to 0.59. In the high d limit, it yields 

This expansion is consistent with the work of Saint James and de Gennes p5fl who got 
the first term of this expansion in the case of an infinite plane in one dimension. In 
two dimensions, one has to take into account the curvature of the cross section. In the 
case of the disc, the equivalent of expansions (||) and (||) have been carried out by Q 
in the limit nd large. 

The hysteresis phenomenon has been rigorously analyzed in |24[] for the supercon- 
ducting solution (n = 0). 

Some interesting questions to be addressed here are: what are the critical fields for 
the vortex solutions losing their local stability (/i* and /i")? can we find them in some 
asymptotic limit such as d large or k large or nd large? and can we prove the existence 
of these hysteresis phenomena for the general case? 



Curve K2(d): 



For fixed k above 0.7, when d is increased from 0, the point (d, k) is first in Region 
1. Then it reaches the critical value dz^n). For d = (^(k), the bifurcation diagram 
1 1 tp\ |oo vs h is decreasing and the superconducting solution bifurcates from the normal 
solution at h = h*. For d a little bigger than ^(k), there is a vortex solution bifurcating 
from the normal solution close to h*. Hence for d = (^(k), at h = /i*, the linearized 
problem near the normal state has two eigenfunctions: one without vortices and one 
with a vortex. This is how uniqueness of solution is lost when increasing d, though the 
assertion needs to be proved mathematically. 

Let D be the fixed domain such that f2 = dD. Then a bifurcated solution near the 
normal state (0, ZioOo) (where ao is such that curl ao = 1 in D and ao • n = on dQ) is 
of the form (ecj), h^a® + eB). Let oo = hod. The second variation of Q near the normal 
state is 



de 2 



Let 



(e<j), h a + eB) = [ \(iV + w 2 a o )0| 2 - n 2 d 2 \(j)\ 2 + |curl B\ 
a Jd 

\{uj) = inf(^|(*V + w 2 a )</>| 2 , [1^2 = 1, G H X {D,C) 



Hence if A(u) > n 2 d 2 , the normal solution is stable, if A(cj) < n 2 d 2 , the normal 
solution is unstable and if X(lo) = K 2 d 2 , degenerate stability occurs. For the eigenvalue 
A(u;) = K 2 d 2 , bifurcation of non normal solutions occurs. Thus, one has to study 

/ (V - iu 2 a ) 2 <j) = X(u)<f> in D, 

1S = on 3D, ^ 

with X(u>) = K 2 d 2 . For most values of k and d, the field h such that the first eigenvalue 
A(cj) is equal to K 2 d 2 yields a single eigenfunction. In region 1, the eigenfunction has 
no vortex while in region 2, the eigenfunction has a vortex. Thus, the curve K2(d) 
corresponds to those values of k and d for which the eigenvalue has two different 
eigenfunctions, one without vortices and one with a vortex. That is, on K,2(d), the 
vortex state starts to exist. 

This situation has been studied in the case of a ball in ||. See also [Q] for more 
recent developments. In this case, the solutions of (Q) with n vortices are of the form 
^ n (r)exp(in9) and have eigenvalue \(u>,n). In particular in Q they draw the function 
\{uj,n) versus uj. The curves X(u>, 0) and \(u>, 1) intersect exactly once for u> = ujq and 
A = Ao- Because of the bifurcation condition X(lo) = K 2 d 2 , it implies that K 2 d 2 = Ao, 
hence the curve K2(d) is of the form nd = constant. It would be interesting to give a 
rigorous proof that the curves A(w, 0) and X(lo, 1) intersect only once for the case of the 
disc and for the case of a more general domain. In region 1, that is below K2(d), the 
first eigenfunction is simple and leads to a solution without vortex. In region 2, that is 
above K2(d), the first eigenfunction is simple and leads to a solution with one vortex, 
but we expect that there is also an eigenfunction with no vortex for a lower field h. 

Similarly, the curves A(u;, n) and X(iv, n + 1) also intersect only once on the numerics 
of H which means in our setting that there are curves nd = C n at which the eigenvalue 
has two eigenfunctions with n and n + 1 vortices. Above nd = C n , a solution with 
n + 1 vortices starts to bifurcate from the normal state and below it, a solution with n 
vortices starts to bifurcate from the normal state, so that the curve nd = C n are the 
critical curves for the existence of n + 1 vortices. 
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In the general case that we are studying, it is totally open to prove that there is 
a unique value of uj such that has two eigenvalues, one with a vortex and one 

without. This would yield to K2(d) = C/d, which is what we have found numerically. 
Moreover, we observe that the field of bifurcation h* satisfies h^d = u hence is constant 
along K2{d). 

Taking this analysis of bifurcation a little further allows us to define 

H(k, d) = {h, s.t. \{Vhd) = K 2 d 2 }. 

In region 1, we expect that Ti has a single element while in Region 2, we expect this 
set to have several elements corresponding to the various branches of solutions with 
several vortices bifurcating from the normal state. But this analysis is open even in 
the case of the disc. 



Region 3: K±{d) < k < Ks(d). 

For parameters in this region, that is large domains and intermediate k (in a relative 
sense), a typical phase diagram is illustrated in Figure 10 with the energy in Figure 
11 and the magnetization in Figure 12. Three solution branches are shown which 



represent the normal solution, the superconducting solution (solid line) and a solution 
with a single vortex (dashed line) . A profile for one of the vortex solutions is given in 
Figure H. 
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Figure 10: The bifurcation curve for d = 6.274 n = 0.35 



The superconducting solution displays a hysteresis behavior as before, but when 
increasing the field, it turns normal instead of going on the vortex branch as in Region 
2. More precisely, as the field is increased, the superconducting solution loses its 
global stability, then its local stability and its drops to the normal branch when the 
super-heating field is reached. Conversely when decreasing the field from the normal 
state, the solution gets on the vortex branch, though it is only a local minimizer. So 
the transition when going down the field is of second order as can be seen on the 
magnetization curve (Figure |l^), but when going up, there is a hysteresis. When 
decreasing the field further, the solution jumps to the superconducting branch. The 
vortex solution is never a global minimizer. In fact, our numerical experiments indicate 
that there are only locally stable vortex solutions. 
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Figure 11: The energy for d = 6.274 k = 0.35 
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Figure 12: The magnetization for d = 6.274 k = 0.35 

We note that, for values of the Ginzburg-Landau parameter k in this region, the 
vortex state has not been frequently studied in the literature, except for superconduct- 
ing sample with extreme geometrical conditions such as thin films, disks and rings. In 
the latter cases, the material displays typical type-II behavior for all ranges of k as the 
Ginzburg-Landau models can be simplified to allow an almost uniform penetration of 
the magnetic field ||. However, the current study is done for three dimensional infinite 
cylinders and the simplified models are not directly applicable. In fact, from the plot 
of the magnetic field given in Figure |l|, we see that there is considerable variation in 
the field strength over the cross section. 

Curve nz{d): 

Let us call H c the thermodynamic critical field introduced by Ginzburg |2(| : the energy 
of the superconducting solution is equal to the energy of the normal solution at this 
field (in our nondimensionalization, it means that the energy of the superconducting 
solution is zero). The curve H3(d) corresponds to the situation where there is a small 
amplitude vortex solution bifurcating from the normal solution exactly at H c . One 
could hope to determine this curve mathematically. 
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We notice that as d tends to infinity, Ks(d) tends to a finite limit close to 0.4 which 
is less than l/v2. Using (||), one can get that in the high d limit, k is close to \\H C , 
where Ai is the first eigenvalue of the linearized problem (0) in an infinite domain and 
is close to 0.59. In Ginzburg's computations (2(J (see also H c is exactly l/\/2 in 

the limit d = oo. This gives Ai/\/2 as a limit for Ks(d) at infinity, which is close to 
0.42. A rigorous mathematical justification of this asymptotic behavior remains to be 
provided. 

Let us recall that the critical value of k equal to l/v2 that separates type I and 
type II superconductors was discovered by Ginzburg studying the bifurcation of the 
superconducting solution from the normal state for an infinite domain. Then Saint 
James and de Gennes discovered surface superconductivity, that is superconductivity 
is nucleated first in the surface and thus appears for higher fields than was calculated 
by Ginzburg. Here, we see that for k less than K^(d), the vortex solution is no longer 
stable but is locally stable near the bifurcation. 



Region 4: k < K±(d) and d > d\{n). 

For parameters in this region, that is k small but domains large enough, the typical 



bifurcation diagram is illustrated in Figure 13 with the energy in Figure 14 and the 



magnetization in Figure 15. 
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Figure 13: The bifurcation curve for d = 4 k = 0.3 



There are superconducting solutions displaying a hysteresis phenomenon and no lo- 
cally stable vortices. The superconducting solution is not always the global minimizer, 
but when increasing the field, the sample remains superconducting until reaching a 
super- heating field h*, where the solution becomes normal with a discontinuous tran- 
sition. Similarly, when decreasing the field, the sample stays normal until the field h* 
which is less than h* , where it turns superconducting by a discontinuous transition. 
Mathematically, we believe that there is a range of fields, between h* and h* where 
there are multiple superconducting solutions. This is in analogy with what happens for 
the one dimensional case. For more rigorous analysis of the one dimensional models, 
we refer to ||]. 

Next, we also note that there is no locally stable vortex solution in the region. It 
is well-known that asymptotically for small k, the vortex solution is not energetically 
favorable, and the material belongs to the typical type-I regime [19, 31 1 where the phase 
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Figure 14: The energy for d = 4 k = 0.3 
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Figure 15: The magnetization for d = 4 k — 0.3 

transition is characterized by the superconducting/normal interface rather than the 
vortex state. When varying the field, the superconducting or normal solutions will not 
turn into a one vortex solution. However, if we do continuation from a vortex solution 
for a bigger value of k (in Region 3), and decrease k, we can still find existence of 
solutions with vortices when we reduce k to values in Region 4 despite of the instability 
of vortex solutions. For even smaller k, continuations in k or in other parameters from 
the vortex solutions fails to produce any new vortex solutions. We believe that when 
decreasing n, the vortex solution first loses its stability near the normal solution (on 
the curve K±(d), but it remains locally stable for Halloo a little higher in the branch. 
For very small k, (especially less than C/d), we believe that there is no vortex solution 
at all. This has been proved in || in the case of a disc. It is an open problem to prove 
that for fixed d and for k small enough, vortex solutions do not exist. 

Curve Ki(d): 

If n is fixed below 0.7 and d is increased, then the curve d\{n) is crossed. It remains 
establish the mathematical existence of this curve. In the particular case where k 
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is very small, the order parameter i/j is almost constant, there are no vortices in the 
domain so that system ([!]) simplifies to 

AA = |V>| 2 A in n curl A = 1 on dfl 

where is a constant that depends on ho. The boundary condition for -0 yields 
/ |^|(|V>| 2 + h^A 2 — 1) = 0, which is the equation of the bifurcation curve. One has 
to find the critical d for which the curve |^|(/io) changes direction of bifurcation near 
| if)\ = 0, so that the bifurcation goes from stable for small d to unstable for larger d. 
Another way to study this curve is to make the bifurcation analysis near the normal 
state, described in the analysis of K2{d) : which yields to (Q). Then one would need 
to take this development to higher order to get the sign of energy of the bifurcated 
branch. This sign changes on K\{d). The fact that the bifurcation from the normal 
state is unstable for large d has not been studied. 

Curve Ki(d): 

Another open problem is to determine the behavior of K^(d) as d tends to infinity. We 
expect it to be of the order C/d for some constant C. We believe that the analysis that 
we have explained for the curve K2(d) is the same here, that is on K^(d) as well the 
eigenvalue has two eigenf unctions. The same analysis of the linearized problem needs 
to be performed. The difference with K%(d) is that the eigenfunction with no vortex is 
stable on K^{d) and unstable on K±{d). 

The point of intersection 

Note that all curves Ki(d) intersect at the same point. Indeed the point of intersection 
of K2(d) and K&{d) has an eigenvalue with two eigenfunctions, one of which (the one 
without vortices) changes stability. Hence this point also belongs to K\(d) since on 
Ki(d) the stability of the solution without vortex changes. Finally this point belongs 
to K^(d) since the energy of the bifurcated solution is zero for both eigenfunctions, in 
particular for the vortex solution. We want to point out that a similar analysis for the 
intersection of these curves has been performed in [Q] in the one dimensional setting 
and it yields to a solvability condition of fourth order at the point of intersection. 

5 Conclusion 

We have obtained very detailed bifurcation diagrams for the Ginzburg-Landau model 
of a two dimensional cross section of a three dimensional superconducting cylinder 
when the applied field is along the direction of the axis. Detailed analysis are provided 
for the solution branches in the different regions of the (k, d) plane. What distinguishes 
the different regions are features like the existence (or the lack of existence) of vortex 
solutions, the global and local stability of solutions, and the hysteresis phenomena. We 
note that the analysis includes regions with small and intermediate values of k which 
is not often featured in the existing studies, as most of the works in the literature focus 
on the vortex state which appears for larger values of k or for thin films. Though we 
have computed the bifurcation diagrams with rectangular cross sections, they are very 
representative of the general cases. Naturally, for small samples, the geometric condi- 
tions may have a stronger effect on leading to detailed alternations to the bifurcation 
curves. A change of topology, such as rings or shells, may present other complications 
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IH but we expect similar bifurcation diagrams remain valid. In addition, comparisons 
with existing theories have been made in the paper. Some remaining questions have 
also been raised. 
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